*****************************************************************;
***figure9.sas                                                ***;
*****************************************************************;
***This program creates data for Figure 9 of the Davis, Grim, ***;
***Haltiwanger, and Streitwieser 2013 RESTAT paper.           ***;
*****************************************************************;
***Output:                                                    ***;
*** ~/csv/figure9_pe.csv                                      ***;
*** ~/csv/figure9_tvs.csv                                     ***;
*****************************************************************;

options ls=80 ps=max spool mprint symbolgen mlogic;

****MACROS*****;

*********************************************************;
***This program is a macro to create dummy variables. ***;
***Required Inputs:                                   ***;
***  ds = dataset name                                ***;
***  dvar = variable for which you want to create     ***;
***         dummies (e.g., industry, state)           ***;
***  prefix = prefix for your dummy (e.g., ind, st)   ***;
***                                                   ***;
***The resulting dummies are numbered from 1 to N     ***;
***(where N is the number of different values)        ***;
***For example, if the prefix is ind, the dummies     ***;
***will be ind1, ind2, ..., indN                      ***;
*********************************************************;
***The dummies are merged back to your original       ***;
***dataset.                                           ***;
*********************************************************;
***WARNING:  Do not run this for a permanent dataset  ***;
***          unless you want the dummies appended to  ***;
***          your permanent dataset.                  ***;
*********************************************************;
***Example:  %makedum(ds=test, dvar=naics, prefix=ind)***;
*********************************************************;

%macro makedum(ds=, dvar=, prefix=);

%*Sort - create dataset with one of each;

proc sort data=&ds nodupkey out=list;
	by &dvar;
run;

%*Count;

data list;
	set list;
	uno = 1;

proc summary data=list;
	output out=count
	sum(uno) = num
	;
run;

%*Put number calculated above into a global macro variable;

data _null_;
	set count;	
	call symput('num', num);
run;

%*Create Dummies;

data list;
	set list (keep = &dvar);
	&prefix.num = _N_;
	%do i=1 %to &num;
		if _N_ = &i then &prefix&&i = 1;
		else &prefix&&i = 0; 
	%end;
run;


%*Merge back to main dataset;

proc sort data=&ds;
	by &dvar;
run;

data &ds;
	merge &ds list;
	by &dvar;
run;

proc datasets lib=work;
	delete list count;
run;
quit;

%mend makedum;

/************************************************************************/
/*                                                                      */
/* Macro: words                                                         */
/* Purpose: determines the number of blank separated words in a         */
/* macro variable or literal.                                           */
/* Author: SAS Guide to Macro Processing, Version 6, 2nd ed, p256       */
/* Usage: %words(&macvar) or %words(literal)                            */
/* Example: %let macvar=This is a test;                                 */
/* %let nwords=%words(&macvar);                                         */
/* %put Number of words in 'This is a test' = &nwords;                  */
/* Number of words in 'This is a test' = 4                              */
/* data _null_;                                                         */
/* nwords = %words(&macvar);                                            */
/* put "Number of words in 'This is a test' = " nwords;                 */
/* Number of words in 'This is a test' = 4                              */
/*                                                                      */
/************************************************************************/

%macro words(string);
%local count word;
%let count=1;
%let word=%qscan(&string,&count,%str( ));
%do %while(&word ne);
%let count=%eval(&count + 1);
%let word=%qscan(&string,&count,%str( ));
%end;
%eval(&count-1)
%mend words;

*Macro to cycle through LHS variables;

%macro lhs(lvar, slvar, wtvar, wtshort);

%r1(pf1, 1, )
%r1(share_gen_hy share_gen_nu share_gen_pn, 2,)
%r1(puc_elect, 3, )
%r1(fac_i_adj, 4, )
%r1(rank_s_adj, 5, )
%r1(pf1 share_gen_hy share_gen_nu share_gen_pn, 6, )
%r1(pf1 share_gen_hy share_gen_nu share_gen_pn puc_elect, 7, )
%r1(pf1 share_gen_hy share_gen_nu share_gen_pn fac_i_adj, 8, )
%r1(pf1 share_gen_hy share_gen_nu share_gen_pn rank_s_adj, 9, )
%r1(pf1 share_gen_hy share_gen_nu share_gen_pn puc_elect rank_s_adj, 10, )
%r1(pf1 share_gen_hy share_gen_nu share_gen_pn puc_elect fac_i_adj rank_s_adj, 11, )
%r1(pf1 share_gen_hy share_gen_nu share_gen_pn puc_elect fac_i_adj int_pf1_fac rank_s_adj, 12, )
%r1(yr1 yr2 yr3 yr4 yr5 yr6 yr7 yr8 yr9 yr10 yr11 yr12 yr13 yr14 yr15 yr16 yr17 yr18 yr19 yr20 yr21
	yr22 yr23 yr24 yr25 yr26 yr27 yr28 yr29 yr30, 13, )
%r1(pf1 share_gen_hy share_gen_nu share_gen_pn puc_elect fac_i_adj int_pf1_fac rank_s_adj 
	yr1 yr2 yr3 yr4 yr5 yr6 yr7 yr8 yr9 yr10 yr11 yr12 yr13 yr14 yr15 yr16 yr17 yr18 yr19 yr20 yr21 yr22 
	yr23 yr24 yr25 yr26 yr27 yr28 yr29 yr30, 14, )
%r1(st1 st2 st3 st4 st5 st6 st7 st8 st9 st10 st11 st12 st13 st14 st15 st16 st17 st18 st19
    st20 st21 st22 st23 st24 st25 st26 st27 st28 st29
    st30 st31 st32 st33 st34 st35 st36 st37 st38 st39 
    st40 st41 st42 st43 st44 st45 st46 st47 st48 st49, 15, )
%r1(pf1 share_gen_hy share_gen_nu share_gen_pn puc_elect fac_i_adj int_pf1_fac rank_s_adj
    st1 st2 st3 st4 st5 st6 st7 st8 st9 st10 st11 st12 st13 st14 st15 st16 st17 st18 st19
    st20 st21 st22 st23 st24 st25 st26 st27 st28 st29
    st30 st31 st32 st33 st34 st35 st36 st37 st38 st39 
    st40 st41 st42 st43 st44 st45 st46 st47 st48 st49, 16, )
%r1(pf1 share_gen_hy share_gen_nu share_gen_pn puc_elect fac_i_adj rank_s_adj 
	yr1 yr2 yr3 yr4 yr5 yr6 yr7 yr8 yr9 yr10 yr11 yr12 yr13 yr14 yr15 yr16 yr17 yr18 yr19 yr20 yr21 yr22 
	yr23 yr24 yr25 yr26 yr27 yr28 yr29 yr30, 17, )
%r1(pf1 share_gen_hy share_gen_nu share_gen_pn puc_elect fac_i_adj rank_s_adj index1_yr_&wtshort, 18, )
%r1(pf1 share_gen_hy share_gen_nu share_gen_pn puc_elect fac_i_adj rank_s_adj index1_yr_&wtshort 
	yr1 yr2 yr3 yr4 yr5 yr6 yr7 yr8 yr9 yr10 yr11 yr12 yr13 yr14 yr15 yr16 yr17 yr18 yr19 yr20 yr21 yr22 
	yr23 yr24 yr25 yr26 yr27 yr28 yr29 yr30, 19, )
%r1(pf1 share_gen_hy share_gen_nu share_gen_pn puc_elect fac_i_adj rank_s_adj index1_yr_&wtshort 
	int_pf1_index1_&wtshort int_hy_index1_&wtshort int_nu_index1_&wtshort int_pn_index1_&wtshort 
	int_elect_index1_&wtshort int_fac_index1_&wtshort int_rank_index1_&wtshort, 20, )
%r1(pf1 share_gen_hy share_gen_nu share_gen_pn puc_elect fac_i_adj rank_s_adj index1_yr_&wtshort 
	int_pf1_index1_&wtshort int_hy_index1_&wtshort int_nu_index1_&wtshort int_pn_index1_&wtshort
        int_elect_index1_&wtshort int_fac_index1_&wtshort int_rank_index1_&wtshort
	yr1 yr2 yr3 yr4 yr5 yr6 yr7 yr8 yr9 yr10 yr11 yr12 yr13 yr14 yr15 yr16 yr17 yr18 yr19 yr20 yr21 yr22 
	yr23 yr24 yr25 yr26 yr27 yr28 yr29 yr30, 21, )

data est_&slvar;
	merge out1 out2 out3 out4 out5 out6 out7 out8 out9 
		out10 out11 out12 out13 out14 out15 out16
		out17 out18 out19 out20 out21
	;
	by vname;

proc export data=est_&slvar outfile="csv/table5_&slvar..csv"
	dbms=csv replace;
	

%*Create dataset with coefficients from all models;

data styr010_c_&slvar;
	set c1 c2 c3 c4 c5 c6 c7 c8 c9 c10
		c11 c12 c13 c14 c15 c16 c17
		c18 c19 c20 c21;

%mend lhs;

*Regression macro;

%macro r1(rvar, mnum, int);

%*Weighted;

proc reg data=r tableout adjrsq mse outest=est&mnum;
	weight &wtvar;
	model &lvar = &rvar &int;
	ods output NObs=bign;
	output out=reg&lvar._&mnum;

proc print data=reg&lvar._&mnum;	
	title1 "CHECK RESID reg&lvar._&mnum";

%**Number of obs;

data bign (keep = model&mnum vname);
	length vname $30.;
	set bign;
	if _N_ = 2;
	rename nobsused = model&mnum;
	vname = 'Nobs';

%**Parameter Estimates;

data parms;
	set est&mnum (drop = &lvar);
	if _TYPE_ = 'PARMS';
	rename _ADJRSQ_ = ADJRSQ;
	rename _EDF_ = EDF;	
	rename _IN_ = IN;
	rename _MSE_ = MSE;
	rename _P_ = P_VALUE;
	rename _RMSE_ = RMSE;
	rename _RSQ_ = RSQ;


%*Continue creation of parameter estimates output dataset;

proc transpose data=parms out=temp1;

data temp1 (keep = vname model&mnum);
	length vname $30.;
	set temp1;
	rename col1 = model&mnum;
	vname = _name_;

%**Standard Errors;

***TEMP;

proc print data=est&mnum;
	title1 "TEMP CHECK est&mnum";

%*Determine if Intercept exists on the est ds;

data _null_;
	%let dsid&mnum = %sysfunc(open(est&mnum));
	%let check&mnum = %sysfunc(varnum(&&dsid&mnum,Intercept));
	%let rc = %sysfunc(close(&&dsid&mnum));

	%if &&check&mnum > 0 %then %do;
		%let rvar2 = Intercept &rvar;
		put "est&mnum has an intercept (&rvar2)";
	%end;
	%else %do;
		%let rvar2 = &rvar;
		put "est&mnum has NO intercept (&rvar2)";
	%end;

run;
quit;

%put "TEST:  check&mnum = &&check&mnum";

%*Count the number of RHS variables;

%let rcount = %words(&rvar2);
%put There are &rcount RHS variables (Model &mnum);


%*Parse RHS variables;

%do ii=1 %to &rcount;
	%let par&ii = %scan(&rvar2, &ii, %str( ));
	%put RHS var &ii = (Model &mnum);
%end;

data se;
	set est&mnum (keep = &rvar2 _TYPE_);
	if _TYPE_ = 'STDERR';
	
	%*Rename RHS vars;
	%do ii=1 %to &rcount;
		rename &&par&ii = &&par&ii.._se;
	%end;

	
%********TEMP;

proc print data=se;
	title1 "TEMP est&mnum - se";


proc transpose data=se out=temp2;

data temp2 (keep = vname model&mnum);
	length vname $30.;
	set temp2;
	rename col1 = model&mnum;
	vname = _name_;

%*Create dataset of coefficients to be merged back to state-year data for estimate creation;

data c&mnum (drop = _TYPE_);
	length model $20.;
	set est&mnum (keep = &rvar2 _TYPE_);
	if _TYPE_ = 'PARMS';
	model = "model&mnum";
	
	%*Rename RHS vars;
	%do ii=1 %to &rcount;
		rename &&par&ii = coeff_&&par&ii;
	%end;

%*Put together;

data out&mnum;
	set temp1 temp2 bign;
	label vname = 'Variable Name';
	label model&mnum = 'Value';

proc sort data=out&mnum;
	by vname;

proc print data=out&mnum;
	title1 "Check out&mnum - &lvar";

proc datasets lib=work;
	delete temp1 temp2 se parms est&mnum;
run;
quit;

%mend r1;

%macro d(slvar, wt, swt);

%*Get coefficients;

data c;
	set styr010_c_&slvar;


%*Get actual values;

data a&slvar;
	set actual (keep = postalst year &slvar._actual &wt);

%*Model 18;

data c18;
	set c;
	if model = 'model18';

%*Look at coefficients relevant to model 18;

proc print data=c18;
	var coeff_intercept coeff_pf1 coeff_share_gen_hy coeff_share_gen_nu coeff_share_gen_pn
		coeff_puc_elect coeff_fac_i_adj coeff_rank_s_adj;
	title1 "Coefficients: Model 18";

%*Create predicted values;

data r18 (keep = postalst year &slvar._:);
	if _N_ = 1 then set c18;
	set r;

	*Elasticity estimate - all variables;
	&slvar._all_18 = coeff_intercept + (coeff_pf1 * pf1) + (coeff_share_gen_hy * share_gen_hy)
				+ (coeff_share_gen_nu * share_gen_nu) + (coeff_share_gen_pn * share_gen_pn)
				+ (coeff_puc_elect * puc_elect) 
				+ (coeff_fac_i_adj * fac_i_adj) + (coeff_rank_s_adj * rank_s_adj)
				+ (coeff_index1_yr_&swt * index1_yr_&swt);
	&slvar._policy_18 = coeff_intercept + (coeff_puc_elect * puc_elect) + (coeff_fac_i_adj * fac_i_adj)
				+ (coeff_rank_s_adj * rank_s_adj);

proc sort data=r18;
	by postalst year;

%*Model 20;

data c20;
	set c;
	if model = 'model20';

%*Look at coefficients relevant to model 20;

proc print data=c20;
	var coeff_intercept coeff_pf1 coeff_share_gen_hy coeff_share_gen_nu coeff_share_gen_pn
		coeff_puc_elect coeff_fac_i_adj coeff_rank_s_adj;
	title1 "Select Coefficients: Model 20";

%*Create predicted values;

data r20 (keep = postalst year &slvar._:);
	if _N_ = 1 then set c20;
	set r;

	*Elasticity estimate - all variables;
	&slvar._all_20 = coeff_intercept + (coeff_pf1 * pf1) + (coeff_share_gen_hy * share_gen_hy)
				+ (coeff_share_gen_nu * share_gen_nu) + (coeff_share_gen_pn * share_gen_pn)
				+ (coeff_puc_elect * puc_elect) 
				+ (coeff_fac_i_adj * fac_i_adj) + (coeff_rank_s_adj * rank_s_adj)
				+ (coeff_index1_yr_&swt * index1_yr_&swt)
				+ (coeff_int_pf1_index1_&swt * int_pf1_index1_&swt)
				+ (coeff_int_hy_index1_&swt * int_hy_index1_&swt)
				+ (coeff_int_nu_index1_&swt * int_nu_index1_&swt)
				+ (coeff_int_pn_index1_&swt * int_pn_index1_&swt)
				+ (coeff_int_elect_index1_&swt * int_elect_index1_&swt)
				+ (coeff_int_fac_index1_&swt * int_fac_index1_&swt)
				+ (coeff_int_rank_index1_&swt * int_rank_index1_&swt)
				;
	&slvar._policy_20 = coeff_intercept + (coeff_puc_elect * puc_elect) + (coeff_fac_i_adj * fac_i_adj)
				+ (coeff_rank_s_adj * rank_s_adj) 
				+ (coeff_int_elect_index1_&swt * int_elect_index1_&swt)
                                + (coeff_int_fac_index1_&swt * int_fac_index1_&swt)
                                + (coeff_int_rank_index1_&swt * int_rank_index1_&swt)
				;

proc sort data=r20;
	by postalst year;


%*Model 13;

data c13;
	set c;
	if model = 'model13';

%*Look at coefficients relevant to model 13;

proc print data=c13;
	var coeff_intercept coeff_yr:;
	title1 "Coefficients: Model 13";

%*Create predicted values;

data r13 (keep = postalst year &slvar._:);
	if _N_ = 1 then set c13;
	set r;

	&slvar._yr_13 = coeff_intercept +
		(coeff_yr1 * yr1) + (coeff_yr2 * yr2) + (coeff_yr3 * yr3) + (coeff_yr4 * yr4) +
		(coeff_yr5 * yr5) + (coeff_yr6 * yr6) + (coeff_yr7 * yr7) + (coeff_yr8 * yr8) +
		(coeff_yr9 * yr9) + (coeff_yr10 * yr10) + (coeff_yr11 * yr11) + (coeff_yr12 * yr12) +
		(coeff_yr13 * yr13) + (coeff_yr14 * yr14) + (coeff_yr15 * yr15) + (coeff_yr16 * yr16) +
		(coeff_yr17 * yr17) + (coeff_yr18 * yr18) + (coeff_yr19 * yr19) + (coeff_yr20 * yr20) +
		(coeff_yr21 * yr21) + (coeff_yr22 * yr22) + (coeff_yr23 * yr23) + (coeff_yr24 * yr24) +
		(coeff_yr25 * yr25) + (coeff_yr26 * yr26) + (coeff_yr27 * yr27) + (coeff_yr28 * yr28) +
		(coeff_yr29 * yr29) + (coeff_yr30 * yr30) 
		;	

proc sort data=r13;
	by postalst year;

%*Model 19;

data c19;
	set c;
	if model = 'model19';

%*Look at coefficients relevant to model 19;

proc print data=c19;
	var coeff_intercept coeff_pf1 coeff_share_gen_hy coeff_share_gen_nu coeff_share_gen_pn
		coeff_puc_elect coeff_fac_i_adj coeff_rank_s_adj coeff_index1_yr_&swt coeff_yr:;
	title1 "Coefficients: Model 19";

%*Create predicted values;

data r19 (keep = postalst year &slvar._:);
	if _N_ = 1 then set c19;
	set r;


	&slvar._all_19 = coeff_intercept + (coeff_pf1 * pf1) + (coeff_share_gen_hy * share_gen_hy)
				+ (coeff_share_gen_nu * share_gen_nu) + (coeff_share_gen_pn * share_gen_pn)
				+ (coeff_puc_elect * puc_elect) + (coeff_fac_i_adj * fac_i_adj)
				+ (coeff_rank_s_adj * rank_s_adj) + (coeff_index1_yr_&swt * index1_yr_&swt) + 
				(coeff_yr1 * yr1) + (coeff_yr2 * yr2) + (coeff_yr3 * yr3) + (coeff_yr4 * yr4) +
				(coeff_yr5 * yr5) + (coeff_yr6 * yr6) + (coeff_yr7 * yr7) + (coeff_yr8 * yr8) +
				(coeff_yr9 * yr9) + (coeff_yr10 * yr10) + (coeff_yr11 * yr11) + (coeff_yr12 * yr12) +
				(coeff_yr13 * yr13) + (coeff_yr14 * yr14) + (coeff_yr15 * yr15) + (coeff_yr16 * yr16) +
				(coeff_yr17 * yr17) + (coeff_yr18 * yr18) + (coeff_yr19 * yr19) + (coeff_yr20 * yr20) +
				(coeff_yr21 * yr21) + (coeff_yr22 * yr22) + (coeff_yr23 * yr23) + (coeff_yr24 * yr24) +
				(coeff_yr25 * yr25) + (coeff_yr26 * yr26) + (coeff_yr27 * yr27) + (coeff_yr28 * yr28) +
				(coeff_yr29 * yr29) + (coeff_yr30 * yr30) 
				;
	&slvar._allxyr_19 = coeff_intercept + (coeff_pf1 * pf1) + (coeff_share_gen_hy * share_gen_hy)
				+ (coeff_share_gen_nu * share_gen_nu) + (coeff_share_gen_pn * share_gen_pn)
				+ (coeff_puc_elect * puc_elect) + (coeff_fac_i_adj * fac_i_adj)
				+ (coeff_rank_s_adj * rank_s_adj) + (coeff_index1_yr_&swt * index1_yr_&swt)
				;
	&slvar._policy_19 = coeff_intercept +
				+ (coeff_puc_elect * puc_elect) + (coeff_fac_i_adj * fac_i_adj)
				+ (coeff_rank_s_adj * rank_s_adj);
	&slvar._yr_19 = coeff_intercept + 
		(coeff_yr1 * yr1) + (coeff_yr2 * yr2) + (coeff_yr3 * yr3) + (coeff_yr4 * yr4) +
		(coeff_yr5 * yr5) + (coeff_yr6 * yr6) + (coeff_yr7 * yr7) + (coeff_yr8 * yr8) +
		(coeff_yr9 * yr9) + (coeff_yr10 * yr10) + (coeff_yr11 * yr11) + (coeff_yr12 * yr12) +
		(coeff_yr13 * yr13) + (coeff_yr14 * yr14) + (coeff_yr15 * yr15) + (coeff_yr16 * yr16) +
		(coeff_yr17 * yr17) + (coeff_yr18 * yr18) + (coeff_yr19 * yr19) + (coeff_yr20 * yr20) +
		(coeff_yr21 * yr21) + (coeff_yr22 * yr22) + (coeff_yr23 * yr23) + (coeff_yr24 * yr24) +
		(coeff_yr25 * yr25) + (coeff_yr26 * yr26) + (coeff_yr27 * yr27) + (coeff_yr28 * yr28) +
		(coeff_yr29 * yr29) + (coeff_yr30 * yr30) 
		;	

proc sort data=r19;
	by postalst year;

%*Model 21;

data c21;
	set c;
	if model = 'model21';

%*Look at coefficients relevant to model 21;

proc print data=c21;
	var coeff_intercept coeff_pf1 coeff_share_gen_hy coeff_share_gen_nu coeff_share_gen_pn
		coeff_puc_elect coeff_fac_i_adj coeff_rank_s_adj;
	title1 "Select Coefficients: Model 21";

%*Create predicted values;

data r21 (keep = postalst year &slvar._:);
	if _N_ = 1 then set c21;
	set r;

	*Elasticity estimate - all variables;
	&slvar._all_21 = coeff_intercept + (coeff_pf1 * pf1) + (coeff_share_gen_hy * share_gen_hy)
				+ (coeff_share_gen_nu * share_gen_nu) + (coeff_share_gen_pn * share_gen_pn)
				+ (coeff_puc_elect * puc_elect) 
				+ (coeff_fac_i_adj * fac_i_adj) + (coeff_rank_s_adj * rank_s_adj)
				+ (coeff_index1_yr_&swt * index1_yr_&swt)
				+ (coeff_int_pf1_index1_&swt * int_pf1_index1_&swt)
				+ (coeff_int_hy_index1_&swt * int_hy_index1_&swt)
				+ (coeff_int_nu_index1_&swt * int_nu_index1_&swt)
				+ (coeff_int_pn_index1_&swt * int_pn_index1_&swt)
				+ (coeff_int_elect_index1_&swt * int_elect_index1_&swt)
				+ (coeff_int_fac_index1_&swt * int_fac_index1_&swt)
				+ (coeff_int_rank_index1_&swt * int_rank_index1_&swt)
		                + (coeff_yr1 * yr1) + (coeff_yr2 * yr2) + (coeff_yr3 * yr3) + (coeff_yr4 * yr4) +
                   		  (coeff_yr5 * yr5) + (coeff_yr6 * yr6) + (coeff_yr7 * yr7) + (coeff_yr8 * yr8) +
               			  (coeff_yr9 * yr9) + (coeff_yr10 * yr10) + (coeff_yr11 * yr11) + 
				  (coeff_yr12 * yr12) + (coeff_yr13 * yr13) + (coeff_yr14 * yr14) + 
				  (coeff_yr15 * yr15) + (coeff_yr16 * yr16) + (coeff_yr17 * yr17) + 
				  (coeff_yr18 * yr18) + (coeff_yr19 * yr19) + (coeff_yr20 * yr20) +
                		  (coeff_yr21 * yr21) + (coeff_yr22 * yr22) + (coeff_yr23 * yr23) + 
				  (coeff_yr24 * yr24) + (coeff_yr25 * yr25) + (coeff_yr26 * yr26) + 
				  (coeff_yr27 * yr27) + (coeff_yr28 * yr28) + (coeff_yr29 * yr29) + 
				  (coeff_yr30 * yr30);
	&slvar._allxyr_21 = coeff_intercept + (coeff_pf1 * pf1) + (coeff_share_gen_hy * share_gen_hy)
				+ (coeff_share_gen_nu * share_gen_nu) + (coeff_share_gen_pn * share_gen_pn)
				+ (coeff_puc_elect * puc_elect) 
				+ (coeff_fac_i_adj * fac_i_adj) + (coeff_rank_s_adj * rank_s_adj)
				+ (coeff_index1_yr_&swt * index1_yr_&swt)
				+ (coeff_int_pf1_index1_&swt * int_pf1_index1_&swt)
				+ (coeff_int_hy_index1_&swt * int_hy_index1_&swt)
				+ (coeff_int_nu_index1_&swt * int_nu_index1_&swt)
				+ (coeff_int_pn_index1_&swt * int_pn_index1_&swt)
				+ (coeff_int_elect_index1_&swt * int_elect_index1_&swt)
				+ (coeff_int_fac_index1_&swt * int_fac_index1_&swt)
				+ (coeff_int_rank_index1_&swt * int_rank_index1_&swt)
				;
	&slvar._policy_21 = coeff_intercept + (coeff_puc_elect * puc_elect) + (coeff_fac_i_adj * fac_i_adj)
				+ (coeff_rank_s_adj * rank_s_adj) 
                                + (coeff_int_elect_index1_&swt * int_elect_index1_&swt)
                                + (coeff_int_fac_index1_&swt * int_fac_index1_&swt)
                                + (coeff_int_rank_index1_&swt * int_rank_index1_&swt)
				;

	&slvar._yr_21 = coeff_intercept + 
		(coeff_yr1 * yr1) + (coeff_yr2 * yr2) + (coeff_yr3 * yr3) + (coeff_yr4 * yr4) +
		(coeff_yr5 * yr5) + (coeff_yr6 * yr6) + (coeff_yr7 * yr7) + (coeff_yr8 * yr8) +
		(coeff_yr9 * yr9) + (coeff_yr10 * yr10) + (coeff_yr11 * yr11) + (coeff_yr12 * yr12) +
		(coeff_yr13 * yr13) + (coeff_yr14 * yr14) + (coeff_yr15 * yr15) + (coeff_yr16 * yr16) +
		(coeff_yr17 * yr17) + (coeff_yr18 * yr18) + (coeff_yr19 * yr19) + (coeff_yr20 * yr20) +
		(coeff_yr21 * yr21) + (coeff_yr22 * yr22) + (coeff_yr23 * yr23) + (coeff_yr24 * yr24) +
		(coeff_yr25 * yr25) + (coeff_yr26 * yr26) + (coeff_yr27 * yr27) + (coeff_yr28 * yr28) +
		(coeff_yr29 * yr29) + (coeff_yr30 * yr30) 
		;	

proc sort data=r21;
	by postalst year;

data rm&slvar;
	merge r18 r13 r19 r20 r21 a&slvar;
	by postalst year;

proc print data=rm&slvar (obs=5);
	title1 "Check rm&slvar";

%*Mean across states by year;

proc sort data=rm&slvar;
	by year;

proc summary data=rm&slvar vardef=wdf;
	by year;
	weight &wt;
	output out=yr&slvar
	mean(&slvar._all_18) = &slvar._all_18
	mean(&slvar._policy_18) = &slvar._policy_18
	mean(&slvar._all_20) = &slvar._all_20
	mean(&slvar._policy_20) = &slvar._policy_20
	mean(&slvar._yr_13) = &slvar._yr_13
	mean(&slvar._yr_19) = &slvar._yr_19
	mean(&slvar._all_19) = &slvar._all_19
	mean(&slvar._allxyr_19) = &slvar._allxyr_19
	mean(&slvar._policy_19) = &slvar._policy_19
	mean(&slvar._yr_19) = &slvar._yr_21
	mean(&slvar._all_19) = &slvar._all_21
	mean(&slvar._allxyr_19) = &slvar._allxyr_21
	mean(&slvar._policy_19) = &slvar._policy_21
	mean(&slvar._actual) = &slvar._actual
	;

%*Calculate mean of annual means across years;

proc summary data=yr&slvar;
	output out=gmean&slvar
	mean(&slvar._actual) = gmean_&slvar._actual
	mean(&slvar._all_18) = gmean_&slvar._all_18
	mean(&slvar._policy_18) = gmean_&slvar._policy_18
	mean(&slvar._all_20) = gmean_&slvar._all_20
	mean(&slvar._policy_20) = gmean_&slvar._policy_20
	mean(&slvar._yr_13) = gmean_&slvar._yr_13
	mean(&slvar._yr_19) = gmean_&slvar._yr_19
	mean(&slvar._all_19) = gmean_&slvar._all_19
	mean(&slvar._allxyr_19) = gmean_&slvar._allxyr_19
	mean(&slvar._policy_19) = gmean_&slvar._policy_19
	mean(&slvar._yr_21) = gmean_&slvar._yr_21
	mean(&slvar._all_21) = gmean_&slvar._all_21
	mean(&slvar._allxyr_21) = gmean_&slvar._allxyr_21
	mean(&slvar._policy_21) = gmean_&slvar._policy_21
	;

data yr&slvar (drop = _TYPE_ _FREQ_);
	if _N_ = 1 then set gmean&slvar;
	set yr&slvar;

	f_&slvar._all_18 = (gmean_&slvar._actual - gmean_&slvar._all_18) + &slvar._all_18;
	f_&slvar._policy_18 = (gmean_&slvar._actual - gmean_&slvar._policy_18) + &slvar._policy_18;
	f_&slvar._all_20 = (gmean_&slvar._actual - gmean_&slvar._all_20) + &slvar._all_20;
	f_&slvar._policy_20 = (gmean_&slvar._actual - gmean_&slvar._policy_20) + &slvar._policy_20;
	f_&slvar._yr_13 = (gmean_&slvar._actual - gmean_&slvar._yr_13) + &slvar._yr_13;
	f_&slvar._yr_19 = (gmean_&slvar._actual - gmean_&slvar._yr_19) + &slvar._yr_19;
	f_&slvar._all_19 = (gmean_&slvar._actual - gmean_&slvar._all_19) + &slvar._all_19;
	f_&slvar._allxyr_19 = (gmean_&slvar._actual - gmean_&slvar._allxyr_19) + &slvar._allxyr_19;
	f_&slvar._policy_19 = (gmean_&slvar._actual - gmean_&slvar._policy_19) + &slvar._policy_19;
	f_&slvar._yr_21 = (gmean_&slvar._actual - gmean_&slvar._yr_21) + &slvar._yr_21;
	f_&slvar._all_21 = (gmean_&slvar._actual - gmean_&slvar._all_21) + &slvar._all_21;
	f_&slvar._allxyr_21 = (gmean_&slvar._actual - gmean_&slvar._allxyr_21) + &slvar._allxyr_21;
	f_&slvar._policy_21 = (gmean_&slvar._actual - gmean_&slvar._policy_21) + &slvar._policy_21;

%*Check;

proc summary data=yr&slvar;
	output out=fmean&slvar
	mean(f_&slvar._all_18) = fmean_&slvar._all_18
	mean(f_&slvar._policy_18) = fmean_&slvar._policy_18
	mean(f_&slvar._all_20) = fmean_&slvar._all_20
	mean(f_&slvar._policy_20) = fmean_&slvar._policy_20
	mean(f_&slvar._yr_13) = fmean_&slvar._yr_13
	mean(f_&slvar._yr_19) = fmean_&slvar._yr_19
	mean(f_&slvar._all_19) = fmean_&slvar._all_19
	mean(f_&slvar._allxyr_19) = fmean_&slvar._allxyr_19
	mean(f_&slvar._policy_19) = fmean_&slvar._policy_19
	mean(f_&slvar._yr_21) = fmean_&slvar._yr_21
	mean(f_&slvar._all_21) = fmean_&slvar._all_21
	mean(f_&slvar._allxyr_21) = fmean_&slvar._allxyr_21
	mean(f_&slvar._policy_21) = fmean_&slvar._policy_21
	;

data gmean&slvar;
	if _N_ = 1 then set gmean&slvar;
	set fmean&slvar;

proc print data=gmean&slvar;
	var fmean_&slvar._policy_18 fmean_&slvar._all_18 fmean_&slvar._policy_20 fmean_&slvar._all_20 
		fmean_&slvar._yr_13 
		fmean_&slvar._all_19 fmean_&slvar._allxyr_19 fmean_&slvar._yr_19 fmean_&slvar._policy_19
		fmean_&slvar._all_21 fmean_&slvar._allxyr_21 fmean_&slvar._yr_21 fmean_&slvar._policy_21
		gmean_&slvar._actual;
	title1 "Check grand means";

%mend d;

%macro i(var);

if &var ne . then good_&var = 1;
else good_&var = 0;

%mend i;

****START PROGRAM****;

*Get data;

data r;
	set elec.styr_elas; *Dataset created in figure7.sas;
	*Drop Nebraska - missing many regulatory variables;
	if postalst ne 'NE';

        *Create interaction variables;
        int_pf1_fac = pf1 * fac_i_adj;
	
	int_pf1_index1 = pf1 * index1_yr;
	int_hy_index1 = share_gen_hy * index1_yr;
	int_nu_index1 = share_gen_nu * index1_yr;
	int_pn_index1 = share_gen_pn * index1_yr;
	int_elect_index1 = puc_elect * index1_yr;
	int_fac_index1 = fac_i_adj * index1_yr;
	int_rank_index1 = rank_s_adj * index1_yr;

	int_pf1_index1_pe = pf1 * index1_yr_pe;
	int_hy_index1_pe = share_gen_hy * index1_yr_pe;
	int_nu_index1_pe = share_gen_nu * index1_yr_pe;
	int_pn_index1_pe = share_gen_pn * index1_yr_pe;
	int_elect_index1_pe = puc_elect * index1_yr_pe;
	int_fac_index1_pe = fac_i_adj * index1_yr_pe;
	int_rank_index1_pe = rank_s_adj * index1_yr_pe;

	int_pf1_index1_tvs = pf1 * index1_yr_tvs;
	int_hy_index1_tvs = share_gen_hy * index1_yr_tvs;
	int_nu_index1_tvs = share_gen_nu * index1_yr_tvs;
	int_pn_index1_tvs = share_gen_pn * index1_yr_tvs;
	int_elect_index1_tvs = puc_elect * index1_yr_tvs;
	int_fac_index1_tvs = fac_i_adj * index1_yr_tvs;
	int_rank_index1_tvs = rank_s_adj * index1_yr_tvs;

*Create year dummies;

%makedum(ds=r, dvar=year, prefix=yr)

*Check year dummies;

proc print data=r (obs=10);
	var year yr:;
	title1 "Check year dummies";

*Create state dummies;

%makedum(ds=r, dvar=postalst, prefix=st)

*Check state dummies;

proc print data=r (obs=10);
	var postalst st:;
	title1 "Check state dummies";

*Output regression dataset for later use;

data styr010_r;
	set r;

*Run regressions;

%lhs(mean_elas_nf_pe, enf_pe, sum_pe, pe)
%lhs(mean_elas_nf_tvs, enf_tvs, sum_tvs, tvs)
%lhs(mean_elas_pe, e_pe, sum_pe, pe)
%lhs(mean_elas_tvs, e_tvs, sum_tvs, tvs)

*Get regression dataset;

data r;
	set styr010_r;

data r2;
	set r (keep = share_gen: puc_elect fac_i_adj rank_s_adj yr: pf1 index1_yr index1_yr_pe index1_yr_tvs);

	%i(share_gen_hy)
	%i(share_gen_nu)
	%i(share_gen_pn)
	%i(puc_elect)
	%i(fac_i_adj)
	%i(rank_s_adj)
	%i(pf1)
	%i(index1_yr)
	%i(index1_yr_pe)
	%i(index1_yr_tvs)

proc freq data=r2;
	tables good:;
	title1 "Check for missings - should be none";

proc means data=r2;
	title1 "Check r2";

*Pull out and rename actual values;

data actual (keep = postalst year fipsst enf_pe_actual e_pe_actual enf_tvs_actual e_tvs_actual sum_pe sum_tvs);
	set r;
	rename mean_elas_nf_pe = enf_pe_actual;
	rename mean_elas_pe = e_pe_actual;
	rename mean_elas_nf_tvs = enf_tvs_actual;
	rename mean_elas_tvs = e_tvs_actual;

proc sort data=actual;
	by postalst year;

%d(enf_pe, sum_pe, pe)
%d(e_pe, sum_pe, pe)
%d(enf_tvs, sum_tvs, tvs)
%d(e_tvs, sum_tvs, tvs)

*Base Model - Figure 9;

data temp_pe;
	set yre_pe (keep = year e_pe_actual f_e_pe_all_18 f_e_pe_policy_18 f_e_pe_yr_19);

proc export data=temp_pe outfile="csv/figure9_pe.csv"
	dbms=csv replace;

data temp_tvs;
	set yre_tvs (keep = year e_tvs_actual f_e_tvs_all_18 f_e_tvs_policy_18 f_e_tvs_yr_19);

proc export data=temp_tvs outfile="csv/figure9_tvs.csv"
	dbms=csv replace;


run;

